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ABSTRACT 

This paper presents the results of conceptual level 
thermal analyses of a High Speed Civil Transport 
(HSCT) wing using p-version finite elements. The work 
was motivated by a thermal analysis of a HSCT wing 
structure which showed the importance of radiation 
heat transfer throughout the structure. The analysis also 
showed that refining a traditional finite element mesh to 
accurately capture the temperature distribution on the 
internal structure led to very large meshes with 
unacceptably long execution times. Further study 
indicated using p-version finite elements might improve 
computation performance for this class of problem. 
Methods for determining internal radiation heat transfer 
were then developed and demonstrated on test problems 
representative of the geometry found in an aircraft wing 
structure. 

This paper presents the results of the application of 
these new methods to the analysis of a high speed 
aircraft wing. Results for both a wing box model as 
well as a full wing model are presented. The reduced 
wing box model allows for a comparison of the 
traditional finite element method with mesh refinement 
(/i-reflnement) to the new p-version finite elements 
while the full wing model demonstrates the 
applicability and efficiency of p-version finite elements 
for large models. 

INTRODUCTION 

This effort began with a study of the thermal analysis of 
a high-speed civilian transport vehicle wing to 
determine the current state of thermal analysis 
capabilities for such a problem. That study, which used 
commercially available thermal analysis software, is 


documented in an earlier paper. 1 The primary difficulty 
in the analysis was accurately modeling the radiation 
heat transfer internal to the wing structure. Initial 
analyses using a coarse mesh and including conduction 
and radiation heat transfer were run successfully; 
however, the large element size produced inaccurate 
temperatures for the internal structure. To improve the 
results, the analysis was repeated using a refined mesh. 
The initial run of the refined model considered heat 
transfer throughout the wing by conduction only, and 
produced a solution with reasonable temperature 
distributions on the internal structure. Internal radiation 
exchange was then added to the model, but the solution 
of the model never successfully completed. 

To get temperature predictions for the internal structure, 
including the effects of internal radiation heat transfer, 
a reduced model was generated. The reduced model 
included a single wing-box section of the wing. This 
model was run both with and without internal radiation 
exchange, and the results showed that internal radiation 
exchange had a significant effect on the temperatures of 
the internal structure. 

Thus a method that could analyze a full wing with 
internal radiation was desired, and the application of 
the p-method appeared promising based on the work in 
reference 1. This paper presents the reults of applying 
the methods developed in reference 1 to the thermal 
analysis of a High Speed Civil Transport (HSCT) wing. 

WING-BOX MODEL 

Before analyzing the full wing, a model of a single 
wing-box was developed. The term wing-box refers to 
the open area enclosed by adjacent ribs and spars and 
the corresponding upper and lower surface sections of a 
wing. The ribs and spars make up the internal structure 
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of the wing with the spars running from the fuselage to 
the wing tip, and the ribs running from the leading edge 
of the wing to the trailing edge of the wing. The ribs 
and spars of the full wing model are shown in Figure 1. 
The wing-box model was developed so that results from 
the new /^-element methods could be compared with the 
results obtained from traditional methods. 



Figure 1: High Speed Civil Transport wing geometry 
including skin, ribs, and spars. 


The wing-box, shown in Figure 2, is 60 inches long 
(fore to aft) by 15 inches wide (port to starboard) by 30 
inches high (thickness of wing). The vertical surfaces 
of the wing-box represent the internal structure (ribs 
and spars) and are adiabatic on their external faces. 

The upper and lower surfaces of the wing-box represent 
the upper and lower surfaces of the wing skin and are 
subjected to convective heating due to the airflow over 
the wing. The convection boundary conditions on the 
upper and lower surface were taken from a mid- wing 
location of the full wing model described in the 
following sections. The methods used to generate the 
convection boundary conditions are discussed in the 
following section as well. Convective heating was 
applied to the upper and lower surfaces of the wing-box 
assuming uniform flow over the surfaces. The 
convection coefficients for both the upper and lower 
surfaces are given in Table 1 and are assumed to remain 
constant throughout the trajectory. Since there is little 
difference in the recovery temperatures for the upper 
and lower surfaces, the same time-dependent values 
were used for both and are listed in Table 2. Heat 
transfer due to natural convection internal to the wing- 
box was ignored. 


z 



Figure 2: Geometry of wing-box model. 


Table 1: Convection coefficients for wing-box model. 


Surface 

Convection Coefficient 
BTU/hr-in 2 -°F 

Upper 

0.06 

Lower 

0.84 


Table 2: Convection recovery temperatures for wing- 
box model. 


Time 

Fluid Temperature 

(hours) 

(°F) 

0.00 

59 

0.20 

334 

4.37 

336 

5.00 

61 


The six surfaces forming the wing-box use a 
corrugated panel construction as shown in Figure 3. 

The panels are made from titanium and the thermal 
properties for this material are listed in Table 3. The 
corrugated panels are modeled with isotropic elements 
with an equivalent thickness of 0.065 inches. This 
thickness represents the thickness of the face sheet 
(0.015 inches) plus the thickness of the core sheet (0.05 
inches). All of the internal faces of the wing-box 
radiate to each other as diffuse-gray surfaces, and the 
upper and lower skin surfaces also radiate to space. All 
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surfaces have an emissivity of 0.85, and the upper 
surface has an absorbed heat flux of 1.5 BTU/hr-in 2 due 
to solar heating. The lower surface has a combined 
solar and earth IR heating of 0.89 BTU/hr-in 2 . 



Table 3: Thermal properties used in wing-box model. 


Material 

Titanium 

Density 

0. 1 6(lb/in 3 ) 

Specific heat 

0.135(BTU/lb-°F) 

thermal conductivity 

0.343(BTU/hr-inch-°F) 


Analysis with p-refinement 

The wing-box was analyzed using the new p-element 
methods with a single element for each of the six 
surfaces of the wing-box. Cases were run with p values 
ranging from 1 to 6. Results were obtained for both the 
radiation sub element and integration methods (see 
reference 1). Each analysis was run using a time step of 
0.1, 0.05, and 0.01 hours with temperature results 
output at every time step. Results for the different time 
steps were similar and the results presented here are 
from the runs using a 0.05 hour time step. Figure 4 
shows temperature contours at time = 0.3 hours for the 
case of p=l using the radiation sub-element method and 
Figure 5 shows the same case with p=2 elements. For 
p=l, the radiation sub-element approach reverts to the 
traditional methods similar :o those employed in 
commercial finite element heat transfer software. 
Comparing the contours for the p= 1 case to the p=2 
case, significant differences can be seen because the 
p=l mesh does not adequately model the side walls of 
the wing-box. For the p=l case, the nodes for the side 
wail elements are located on upper and lower surfaces 
with no nodes in-between. With no nodes between the 
upper and lower surfaces, the temperature of the 
sidewalls cannot be computed accurately. The nodes 
on the upper and lower surfaces are shared by the upper 
and lower surface elements which have the convection 


boundary conditions. The convection boundary 
conditions drive the temperatures of these nodes and the 
temperature lag midway between the upper and lower 
surfaces is not captured. 

Figure 6 shows temperature contours for the integration 
method case with p= 4 elements, at 0.3 hours after 
takeoff. The results are similar to the p=2 results in 
Figure 5, with the temperatures in the middle of the 
vertical surfaces (the internal structure) slightly 
warmer. In fact, the results do not change significantly 
for higher-order elements. Figure 7 shows the results 
for the integration method using p^6 elements. In the 
cases run here, the results for the radiation sub element 
method and the integration method were similar. 

To see the transient nature of the temperatures gradients 
in the wing-box, two temperatures from the forward 
rib-port spar junction have been plotted versus time. 
Point A is at the intersection of the fore spar, the port 
rib, and the lower surface while point B is along the 
spar-rib edge halfway between the upper and lower 
surfaces (see Figure 7). Figure 8 shows the temperature 
histories computed using the radiation sub-element 
method with p=2 elements, and Figure 9 shows the 
results from the integration method using the p=6 
results. Figure 8 shows that point B initially drops in 
temperature with the p=2 elements. This behavior is a 
consequence of using a consistent mass matrix in the 
transient solution algorithm in a problem with rapid 
temperature change (see, for example, p. 334 of [2]). 
The temperature drop had the same approximate 
magnitude and duration regardless of the time step used 
(0.1 hours, 0.05 hours, 0.01 hours). Switching to a 
backwards difference for the time derivative eliminated 
the temperature drop; however, the backwards 
difference is only first-order accurate whereas the 
central difference used in the Crank-Nicholson method 
is second-order accurate. Note however, that the 
temperature drop did not occur for the higher-order 
elements as shown Figure 9 for the p^ elements. Lobo 
and Emery [3] studied this behavior and concluded that 
the higher-order elements are not immune from this 
anomalous behavior but simply less susceptible to it. 
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Figure 4: Temperature contours of the wing-box model 
using the radiation sub element method with p=l 
elements at 0.3 hours. 


Figure 6: Temperature contours of the wing-box model 
using the integration method with p=4 elements at 0.3 
hours. 



Figure 5: Temperature contours of the wing-box model 
using the radiation sub-element method with p=2 
elements at 0.3 hours. 


Figure 7: Temperature contours of the wing-box model 
using the integration method with p=6 elements at 0.3 
hours. 
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Figure 8: Transient temperature response at selected 
points of the wing-box model using the radiation sub- 
element method with p~2 elements. (Location of points 
shown in Figure 7.) 


The traditional method does however have one 
advantage, speed. Table 4 lists the number of degrees 
of freedom and the execution time for several cases. 

The execution times were obtained on a 133 MHz 
Pentium PC with 32 Mbytes of RAM. For a given 
number of degrees of freedom the traditional method is 
faster than the radiation sub-element method which in 
turn is faster than the integration method. The 
differences in run times increase significantly as the 
element order increases in the p-method. However, a 
comparison based on the degrees of freedom does not 
account for the improved accuracy of the higher-order 
methods. The contour plots show the p=2 solution to 
be the most similar to the traditional method with a 10 
by 1 0 mesh (the most refined mesh used) but with a 
significantly lower execution time. This result is 
consistent with the results in reference [ 1 ] where a more 
rigorous accuracy comparison of the methods was 
presented. 



Figure 9: Transient temperature response at selected 
points of the wing-box model using the integration 
method with p=6 elements. (Location of points shown 
in Figure 7.) 

Analysis with ^-refinement 

Several finite element analyses were made of the wing- 
box using the traditional radiation methods with h- 
reflnement. In these analyses, each surface was 
subdivided into an n by n mesh where n varied from 1 
to 10 (the view factor code was the limiting factor in 
increasing the mesh size). Results for the 2 by 2, 5 by 
5, and 10 by 10 meshes are shown in Figure 10-12. The 
/j-refinement behavior shown in these figures is similar 
to the behavior seen in the p-enrichment cases, but note 
that the lowest order p results (p=2) are similar to the 
results from the largest n by n mesh case. 



Figure 10: Temperature contours for wing-box model 
case 12 (2 elements per edge or 4 elements per surface) 
at 0.3 hours. 
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Table 4: Execution times for wing-box model using 
new methods versus traditional methods. 



Figure 1 1 : Temperature contours for wing-box model 
case 14 (5 elements per edge, 25 elements per surface) 
at 0.3 hours. 



Figure 12: Temperature contours for wing-box model 
case 15(10 elements per edge, 100 elements per 
surface) at 0.3 hours. 


Case 

Method* 

element 

P 

elements 
per edge 

(n) 

#of 

DOF 

Execu- 
tion time 
(sec.) 

1 

RSE 

2 

1 

26 

72 

2 

RSE 

3 

1 

56 

193 

3 

RSE 

4 

1 

98 

502 


RSE 

5 

1 

152 

1261 

5 

RSE 

6 

1 

218 

2733 

6 

IM 

2 

1 

26 

315 

7 

IM 

3 

1 

56 

996 

8 

IM 

4 

1 

98 

2867 

9 

IM 

5 

1 

152 

6777 

10 

IM 

6 

1 

218 

16200 

11 

TM 

1 

1 

8 

9.1 

12 

TM 

1 

^ 2 

26 

37 

13 

TM 

1 

3 

56 

84 

14 

TM 

1 

5 

152 

360 

15 

TM 

1 

10 

602 

17018 


*RSE = Radiation sub-element method, IM = 
Integration method, TM = Traditional Method 


FULL WING ANALYSIS 

Description 

The next analysis is of a complete HSCT wing. The 
analysis presented here is typical of a conceptual level 
design effort that integrates several disciplines to 
evaluate the feasibility of a proposed vehicle. The 
system analysis begins with an aerodynamic analysis of 
the proposed vehicle at selected points throughout the 
trajectory. Aerodynamic heating data is then used in a 
thermal analysis of the vehicle to predict temperatures 
throughout the structure. These temperatures are then 
transferred to a structural model of the vehicle along 
with the pressure loads from the aerodynamic analysis. 
The structural model is used to size the load-carrying 
components of the vehicle, typically through an 
optimization process. The final component sizing 
determines component level and overall vehicle weight. 
The weight data, trajectory data, and aerodynamic data 
are then used in a trajectory analysis to determine the 
overall viability of the vehicle. The entire process is 
then iterated until a final optimized design is 
determined. 

The conceptual level design used in this process 
produces thermal and structural models with simple 
geometric elements. For example, holes that 
commonly occur in the internal structure of the wing 
are not considered at this stage, but are left to the 
detailed design process. The goal of the thermal and 
structural analyses in this process is not to analyze a 
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particular design in great detail, but rather to accurately 
evaluate multiple conceptual designs in order to arrive 
at an optimum structure. 

The thermal analysis of an entire vehicle is often 
performed using individual models for the major 
components of the vehicle. The analysis presented here 
is of a wing structure of a high-speed civil transport. 
The geometry used in the analysis was supplied by the 
GEOLAB group at NASA Langley Research Center. 
The wing, shown in Figure 1, is approximately 1 13 feet 
long at the root, and 55 feet wide at the trailing edge. 
The structure consists of the upper and lower skin 
surfaces of the wing as well as the internal structure of 
ribs and spars. The skin surfaces represent the outer 
mold line of the wing and are consistent with the 
geometry used for aerodynamic analysis. 

The wing skin uses the same corrugated panel 
construction used in the wing-box model, which is 
shown in Figure 3. The panels were modeled using the 
same material properties used for the wing-box model 
which are listed in Table 3. As in the wing-box model, 
the corrugated panels are modeled with isotropic 
elements with an equivalent thickness of 0.065 inches. 
This thickness represents the thickness of the face sheet 
(0.015 inches) plus the thickness of the core sheet (0.05 
inches). For simplicity, the internal structure was 
assumed to be made of the same corrugated panel 
construction. 

Aerodynamic Heating 

A five-hour flight trajectory representative of a 
commercial airline or transport route was used in the 
analysis. The Mach 2.4 cruise lasts 3.8 hours and 
reaches a maximum altitude of 70,000 feet. Figure 13 
shows the altitude and speed of the vehicle throughout 
the trajectory, the analysis does not include any special 
maneuvers (for example an emergency descent) which 
often drive the design of an aircraft. 





Figure 13: Mach number and altitude for commercial 
route of a High Speed Civil Transport. 


To compute the convection boundary conditions over 
the full wing, NASA Langley’s version of the 
MINIVER computer code, LANMIN [4], was used. 

The upper and lower surfaces were subdivided to give 
one surface element between adjacent pairs of ribs and 
spars. The surface elements of Figure 1 represent the 
surface elements used to compute the convection 
boundary conditions. The centroid of each surface was 
computed along with the distance from the leading edge 
of the wing to the surface centroid. This data was then 
combined with the trajectory data to compute a 
convection coefficient and recovery temperature for 
every surface element in Figure 1 at every trajectory 
point. Upon inspection of the data, and based on results 
of some wing-box model runs, it was determined that 
the convection boundary condition could be accurately 
modeled using only four trajectory points and 
interpolating the data between points. This greatly 
reduced the amount of data that had to be transferred to 
the finite element model without a significant loss in 
accuracy. 

Full Wing Finite Element Thermal Model 

The full wing model was meshed using one element to 
represent each surface shown in Figure 1 . Thus the 
mesh used in the aerodynamic heating analysis is 
coincident with the thermal analysis mesh and each 
wing surface element has a separate transient 
convection boundary condition associated with it. The 
wing was assumed to be totally isolated from the 
fuselage in this model, thus an adiabatic boundary 
condition was assumed at the fuselage. The internal 
surfaces of the wing were assumed to be diffuse-gray 
surfaces with an emissivity of 0.85. Results from the 
wing-box model showed that the convection boundary 
condition dominated the surface heating conditions 
such that external radiation and solar and earth IR heat 
fluxes had negligible effect, so they were not included 
in the full wing analysis. Natural convection inside the 
wing was also not considered. 

Due to the large number of elements, and the problems 
encountered with view factor computations for such 
large models in the initial study, only the integration 
method was used to analyze the full wmg model. The 
model consists of 235 elements and 725 nodes, and was 
analyzed using p=2, p=3, and p=4 elements. There are 
725 degrees of freedom in the model using p=2 
elements, 1769 degrees of freedom in the p=3 element 
model, and 3283 degrees of freedom in the p=4 element 
model. The model with p=2 elements was run with 
time steps of 0.1, 0.05, and 0.01 hours with no 
significant differences in the results. The p=3 and p=4 
runs used a time step of 0.05 hours. All runs were 
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performed on a 200 MHz Pentium class personal 
computer. 

A temperature contour plot of the upper surface of the 
wing at the beginning of the cruise portion of the 
trajectory (0.2 hours into the flight) is shown in Figure 
14. First note that the temperature variation over the 
entire wing surface is not that large, 10 to 15 °F. 

Except for the outboard leading edge, the lowest 
temperatures on the wing occur near the fuselage where 
the wing is thickest. The deeper ribs and spars in this 
area act as heat sinks preventing the surrounding 
structure from heating as fast as other areas. Moving 
out from the fuselage (+y) the wing thins and the 
resulting internal structure becomes smaller. The 
smaller internal structure provides less thermal mass 
and the corresponding skin temperatures increase. Near 
the fuselage, the heat sink effect of the relatively 
massive internal structure can be clearly seen by the 
oval contours corresponding to individual wing-boxes. 
The oval contours occur because the temperature at the 
edges of the wing-boxes are held down by the adjacent 
internal structure while the temperatures at the center of 
the wing-box respond more rapidly. Figure 15 shows 
the corresponding temperature contours for the lower 
surface. Since the lower surface has a higher 
convection coefficient, it heats up faster than the upper 
surface. This is evidenced by the overall higher 
temperature in Figure 15 as compared to Figure 14. 

The temperature gradients on both surfaces disappear 
quickly, and are in fact gone by 0.3 hours as shown in 
Figure 16 and Figure 17. 

To see the temperature response of the internal 
structure, the data for a spar, a rib, and a single wing- 
box have been extracted from the full wing data set. 

The approximate locations of the components are 
shown in Figure 17. Figure 18 shows the temperature 
contours for the selected rib at 0.2 hours. The same 
general patterns that were seen in the wing-box model 
appear here. The upper and lower surfaces of the spar 
are at uniform temperatures while the temperatures 
between the upper and lower surfaces are cooler. The 
lowest temperatures occur where the spars intersects the 
rib (vertical lines in Figure 18) halfway between the 
upper and lower surfaces. These patterns are also 
shown in Figure 19 where the temperatures along the 
upper and lower surface of the spar are plotted along 
with the temperature halfway between the upper and 
lower surfaces. The symbols on the line plot are used 
only to differentiate the curves, they do not represent 
the actual data points (each curve was generated using 
approximately 60 equally spaced temperature values). 

The mid-plane temperature profile shows a steep drop 
at the intersection of each spar. Figure 20 shows the 
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same temperature plots at 0.3 hours. The temperature 
variation along the mid-plane at this later time shows 
the same patterns that were present at 0.2 hours, 
although the overall temperature is higher and the 
magnitude of the variations has dimmished. Figure 2 1 
shows the results from the model with p=4 elements at 
0.2 hours. The results are similar to the p=2 results with 
the p=4 elements producing a flatter temperature profile 
near the center of the wing-boxes. 

The temperature contours for the selected spar at 0.2 
hours using p=2 elements are shown in Figure 22. 

Once again the upper and lower surfaces are relatively 
warm with uniform temperatures. The lowest 
temperatures occur midway between the upper and 
lower surfaces where the ribs intersect the spar. As the 
wing thins moving out from the fuselage, the mid-plane 
temperatures (halfway between the upper and lower 
surfaces) increase and the gradients in the spar 
diminish. This pattern is also illustrated in Figure 23 
where the temperatures along the upper, lower, and 
mid-plane surfaces for this spar at 0.2 hours are plotted. 
Figure 24 shows the same temperature plots at 0.3 
hours. Once again, the overall structure temperature 
has increased and the magnitude of the temperature 
gradients have decreased relative to the earlier time 
point. Figure 25 shows the results from the model with 
p=4 elements at 0.2 hours. Again the results are similar 
to the p=2 results with the p=4 elements producing a 
flatter temperature profile near the center of the wing- 
boxes. 

Finally, Figure 26 shows the temperature contours of 
the selected wing-box using the p=2 results at 0.3 
hours. This wing-box is one of the larger wing-boxes 
in the wing and thus has some of the largest 
temperature gradients. The contours shown in Figure 
26 have the same patterns as seen in the wing-box 
analysis. The coolest temperatures occur at the 
intersection of the ribs and spars halfway between the 
upper and lower surface. One slight difference is that 
the outboard rib is shorter than the inboard rib and 
results in lower gradients across the outboard section of 
the wing-box. 
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Figure 14: Temperature contours of the full wing model Figure 17: Temperature contours of the full wing model 
upper surface at 0.2 hours using p=2 elements. lower surface at 0.3 hours using p=2 elements. 



Note: z dimension scale expanded to 5 times x dimension scale 


Figure 15:Temperature contours of the full wing model Figure 18: Temperature contour for selected rib of full 
lower surface at 0.2 hours using p=2 elements. wing model at 0.2 hours using p=2 elements. 
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Figure 19: Temperatures along selected rib of full wing 
model at 0.2 hours using p=2 elements. 
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Figure 20: Temperatures along selected rib of full wing Figure 23: Temperatures along selected spar of full 

model at 0.3 hours using p=2 elements. wing model at 0.2 hours using p=2 elements. 



1800 2000 2200 2400 _ 

* pncha*] r igure 24: Temperatures along selected spar of full 


Figure 21 : Temperatures along selected rib of full wing W * n ® moc * e * at ^ hours us * n & P 2 elements, 

model at 0.2 hours using p=4 elements. 



Note: z dimension seal* expanded to 5 times x dimension scale 

Figure 22: Temperature contour for selected spar of full 
wing model at 0.2 hours using p=2 elements. 
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Figure 25: Temperatures along selected spar of full 
wing model at 0.2 hours using p=4 elements. 


stage of a vehicle when overall design parameters are of 
interest. The new methods will permit a better heat 
transfer representation of the internal structure and will 
thus give a more accurate prediction of the temperature 
distribution throughout the wing. The improved 
temperature distribution can then be used in a structural 
model of the wing to evaluate thermal stress. The 
refined temperature distribution will improve the 
thermal stress distribution producing a more accurate 
sizing of the structural components of the wing, an 
important factor since these vehicles are extremely 
weight critical. 



Figure 26: Typical wing-box extracted from full wing 
model, time = 0.3 hours, p=2 elements. 
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CONCLUSIONS 

The methods developed in [1] have been used 
for the thermal analysis of a High Speed Civil 
Transport wing using a wing-box model as well as a 
full wing model. Results from the full wing model 
indicate that a solution with p=2 elements provides a 
good balance between accuracy and cost (execution 
time). The development of these methods was the 
result of an initial study that demonstrated the 
difficulties in analyzing a full wing model due to the 
complexities involved with the internal radiation heat 
exchange. The new methods provide a way to analyze 
the full wing structure accurately and efficiently. This 
capability is particularly useful in the conceptual design 
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